This invention relates to geophysical exploration for petroleum and minerals. More particularly, this invention is directed to geophysical prospecting by means of the seismic technique.
Seismic prospecting involves generating seismic waves at the surface of the earth by means of a seismic source. The seismic waves travel downward into the earth and are reflected and/or refracted due to differences in acoustic impedance at the interfaces of various subsurface geological formations. Detectors, called seismometers, or geophones, located along the surface of the earth and/or in a borehole produce analog electrical seismic-trace signals in response to detected seismic wave reflections and/or refractions. The analog electrical seismic-trace signals from the seismometers, or geophones, can then be recorded. Alternatively, the analog electrical seismic-trace signals from the seismometers, or geophones, can be sampled and digitized prior to being recorded. The seismic-trace data recorded in either manner is subsequently processed and analyzed for determining the nature and structure of the subsurface formations. Specifically, this invention is directed to the suppression of noise which is present in the seismic-trace data, especially in the case where a low energy surface seismic source, such as a vibrator, is used for imparting seismic energy to the earth.
Many techniques for generating and recording seismic waves are currently in use. Exploding-gas and compressed-air guns placed on the surface of the earth and dynamite are examples of high energy seismic sources which generate a sharp pulse (impulse) of seismic energy. Vibrators, which generate a "chirp" signal of seismic energy, and hammers are examples of low energy surface seismic sources. In the case of vibrators, the recorded seismic wave reflections and/or refractions are cross-correlated with a replica (called the "pilot signal") of the original "chirp" signal in order to produce recordings similar to those which would have been produced with a high energy impulsive seismic source. This process is known as "vibroseis."
Considered in more detail, vibroseis seismic prospecting, commercialized by Continental Oil Company, typically employs a large, vehicle-mounted vibrator as a seismic source. The vehicle is deployed to a prospect area, and the vibrator is positioned in contact with the surface of the earth. Thereafter, the vibrator is activated for imparting vibrations to the earth, thereby causing seismic waves to propagate through the subsurface formations. The seismic wave reflections and/or refractions are detected by seismometers, or geophones, deployed in the prospect area.
Advantageously, the use of a vibrator can be more economical than the use of dynamite. Furthermore, as compared to the use of a high energy impulsive seismic source, such as dynamite, the frequency of the seismic waves generated by a vibrator can be selected by controlling the frequency of the pilot signal to the power source, such as a hydraulic motor, which drives the vibrator. More particularly, the frequency of the pilot signal to the vibrator power source can be varied, that is, "swept," for obtaining seismic-trace data at different frequencies. Consider, for example, Doty et al. U.S. Pat. No. 2,688,124 which discloses how a low energy seismic wave, such as generated by a vibrator, can be used effectively for seismic prospecting if the frequency of the vibrator "chirp" signal which generates the seismic wave is swept according to a known pilot signal and the detected seismic wave reflections and/or refractions are cross-correlated with the pilot signal in order to produce seismic-trace recordings similar to those which would have been produced with a high energy impulsive seismic source. Tyically, the pilot signal is a swept frequency sine wave which causes the vibrator power source to drive the vibrator for coupling a swept sine wave "chirp" signal into the earth. A typical swept frequency operation can employ, for example, a 10- to 20-second long sine wave "chirp" signal with a frequency sweep of 14 to 56 Hz. The swept frequency operation yields seismic-trace data which enables the different earth responses to be analyzed, thereby providing a basis on which to define the structure, such as the depth and thickness, of the subsurface formations.
Unfortunately, recorded seismic-trace data always includes some background (ambient) noise in addition to the detected seismic waves reflected and/or refracted from the subsurface formations (referred to as the "true earth response"). Noise can be classified as "stationary" and "nonstationary." In probabilistic terms, both stationary and nonstationary noise can be random. Stationary noise is random noise, such as atmospheric electromagnetic disturbances. Stationary noise is statistically time-invariant over the period of acquisition of seismic-trace data for producing a recording. Nonstationary noise is random and often occurs as bursts or spikes generally caused by wind, traffic, recorder electrical noise, etc. Nonstationary noise is statistically time-variant over the period of acquisition of seismic-trace data for producing a recording and exhibits relatively large excursions in amplitude.
When high energy impulsive seismic sources are used, the level of the detected true earth response seismic signal is usually greater than the ambient noise. However, when low energy surface seismic sources are used, the ambient noise can be at a level greater than the true earth response seismic signal. For this reason, seismic-trace recordings are often made involving the repeated initiation of a low energy surface seismic source at about the same origination point, thereby producing a sequence of seismic-trace data based on seismic wave reflections and/or refractions that have traveled over about the same path and therefore have approximately the same travel times. The process of adding such seismic-trace data together for improving the signal-to-noise ratio of the composite seismic-trace recording is known as "vertical compositing" or "vertical stacking." It should be distinguished from "horizontal stacking," a process applied to a sequence of seismic-trace data based on seismic wave reflections from approximately the same subsurface point (referred to as the "common-depth point," or "CDP") but which has been generated and recorded at different surface locations. Horizontal stacking of CDP seismic-trace data requires that time corrections (called "normal moveout," or "NMO," corrections) be applied before the traces are summed together, since travel times from seismic source to detector are unequal for each trace in the sequence. While this invention could be applied in either process, it is primarily intended to improve the vertical stacking process.
In connection with the earlier mentioned swept frequency operation of vibroseis seismic prospecting, it is common practice to vertically stack, or sum, the seismic-trace data from a series of initiations, that is, sequential swept frequency operations, to produce a composite seismic-trace recording for the purpose of improving the signal-to-noise ratio of the seismic-trace data. However, since the earliest days of vibroseis seismic prospecting, which is most economical when conducted along the existing road network where a large vehicle-mounted vibrator can be used, noise, and in particular, nonstationary noise such as burst noise associated with road traffic, has been recognized to have a severe adverse effect on seismic-trace data quality. Unless the nonstationary noise is somehow suppressed vis-a-vis the true earth response seismic signal, the ability to accurately map the subsurface formations is diminished.
Unfortunately, the commonly used technique described above of vertically stacking the seismic-trace data for the purpose of improving the signal-to-noise ratio has proven inadequate in the presence of nonstationary noise which appears during vibroseis seismic prospecting. That is, the low instantaneous transduced signal level of detected seismic wave reflections and/or refractions in the case of vibroseis seismic prospecting requires there be a long vibrator "chirp" signal duration, either a single very long swept frequency "chirp" signal or, more likely, a shorter swept frequency "chirp" signal (10-20 seconds) repeated many times. However, during the swept frequency operation, a large burst of noise will swamp the low instantaneous transduced signal level of detected seismic wave reflections and/or refractions and if digitized and vertically stacked will render the seismic-trace data unusable. The longer the duration of the greater the number of repetitions of the swept frequency "chirp" signal, the greater the risk of exposure to such fatal bursts of noise.
Interestingly, early analog field recording had such limited dynamic range that the noise bursts saturated the recording medium, whereby the noise was moderated to the extent that the recording was not rendered unusable. Digital field recording, on the other hand, with a cableless seismic digital recording system, such as the one disclosed in Weinstein et al. U.S. Pat. No. 3,946,357, records such noise bursts faithfully, thereby rendering the recording unusable. There is a need to improve the signal-to-noise ratio of seismic-trace data collected by digital field recording during seismic prospecting with a low energy surface seismic source in a noisy environment, particularly where nonstationary noise appears.
The following provides a more detailed analysis of known approaches which involve vertical stacking for improving the signal-to-noise ratio of seismic-trace data. Let the j-th digitized sample of the i-th seismic-trace signal (X.sub.i,j) in a sequence which is to be vertically stacked be represented by: EQU X.sub.i,j =.alpha..sub.i (.sigma..sub.i,j +.eta..sub.i,j) (1)
where .sigma..sub.i,j is the true earth response seismic signal, .eta..sub.i,j is the noise, and .alpha..sub.i is a scale factor (scalar) corresponding to seismic source and/or detector earth coupling and recorder amplifier gain variations.
The assumptions can be made that the true earth response seismic signal embedded in each trace is coherent and in phase (correlated) and that the noise is random and incoherent (uncorrelated) with zero mean value. Under such assumptions, the square value (variance) of any N-sample time "window," or portion, of the i-th trace is: ##EQU1##
Since the noise is uncorrelated with the signal, the expected value of the middle term in Equation (2) is zero. In practice, the value is small and therefore can be neglected, resulting in: ##EQU2## where S.sub.i.sup.2 and n.sub.i.sup.2 are the received earth response seismic signal and noise variances, respectively, and .gamma..sub.i is the signal-to-noise power ratio of the i-th trace.
In general, the objective of vertical stacking is to maximize the signal-to-noise ratio of the resultant recording. To accomplish that objective, each seismic-trace signal sample is multiplied by a weighting function or scalar and summed with the other similarly weighted samples in the M-trace sequence. The j-th sample of the vertically stacked recording is then represented by: ##EQU3## where .beta..sub.i is the weighting function or scalar. The signal-to-noise power ratio of the vertically stacked recording is: ##EQU4##
Since the noise is uncorrelated from trace to trace within the M-trace sequence, the expected value of the cross terms in the denominator of Equation (5) is zero. In practice, the value of the summed cross terms is small and therefore can be neglected, resulting in: ##EQU5##
In order to determine the optimum weighting values which maximize the signal-to-noise power ratio, the partial derivative of .gamma. in Equation (6) with respect to some .beta..sub.K, where K is between 1 and M, inclusive, is equated to zero. The result of this operation after some simplification is: ##EQU6##
In the case M=1, Equation (7) is trivial, and .beta..sub.1 is arbitrary.
In the case M=2, we let K=2, and after algebraic mainpulation, we then have EQU .beta..sub.2 =.beta..sub.1 (n.sub.1.sup.2 s.sub.2 /n.sub.2.sup.2 s.sub.1) (8)
Now, we let .beta..sub.1, an arbitrary quantity, be .beta..sub.1 =s.sub.1 /n.sub.1.sup.2. Therefore, Equation (8) becomes: EQU .beta..sub.2 =s.sub.2 /n.sub.2.sup.2 ( 9)
Using mathematical induction, we assume that .beta..sub.i =s.sub.i /n.sub.i.sup.2 for all i between 1 and M-1, inclusive. In order to prove .beta..sub.M =s.sub.M /n.sub.M.sup.2, then it follows from Equation (7), letting K=M: ##EQU7##
After algebraic manipulation, we have: EQU .beta..sub.M =s.sub.M /n.sub.M.sup.2 ( 11)
Consequently, the optimum weighting value for the i-th trace is: EQU .beta..sub.i *=s.sub.i /n.sub.i.sup.2 =.gamma..sub.i /s.sub.i ( 12)
In order to maximize .gamma., Equation (12) requires that each seismic-trace signal sample be weighted in proportion to its true earth response seismic signal amplitude and inversely to its noise power. Substituting this requirement into Equation (6) yields: ##EQU8## which implies that, under optimum weighting, the signal-to-noise power ratio of the optimum vertically stacked recording is equal to the sum of the signal-to-noise power ratios of the seismic-trace signals.
While Equation (12) is mathematically exact under the assumptions of coherent, in-phase true earth response seismic signal and random noise, computation of the optimum weighting value requires statistical estimation of s.sub.i and n.sub.i.sup.2, or s.sub.i and .gamma..sub.i. Among others, Robinson, "Statistically Optimal Stacking of Seismic Data," Geophysics, June 1970, proposes and evaluates schemes for determining these statistical estimates through auto-correlations and cross-correlations of the seismic-trace signals. Such schemes require considerable computer execution time and memory storage, thereby rendering them impractical and uneconomical in field recording environments.
However, Robinson's application to synthetic and actual field seismic-trace data demonstrates that while a maximum signal-to-noise ratio of the vertically stacked seismic-trace data is achieved when statistical estimates of s.sub.i and n.sub.i.sup.2 are available, sufficiently improved results are possible with approximations. The approximations relate only to the manner in which the weighting values are determined.
The simplest approximation occurs when s.sub.i and n.sub.i.sup.2 do not change from 1 seismic-trace signal to the next. In this case, .beta..sub.i in Equation (12) is a constant, and computation of relative weighting values is not required. If s.sub.i and n.sub.i.sup.2 are constant, then so too is .gamma..sub.i. Therefore, the signal-to-noise power ratio improvement in the vertically stacked recording is simply .gamma.*=M.gamma. where .gamma. is the signal-to-noise power ratio of the seismic-trace signals. Also, note that the signal-to-noise amplitude ratio is improved by a factor of .sqroot.M. This approximation, often referred to as "true amplitude" summation, is implemented in various commercially available digital field recorders.
However, consider several repetitions of a swept frequency "chirp" signal at each of a plurality of vibration points for providing a set of seismic-trace data from which the true earth response is to be statistically estimated. Assume one recording within the set to be in the presence of very high burst noise. Because of the nonstationarity of burst noise, true amplitude summation in fact gives an estimate of the noise, not the true earth response, and therefore, a simple addition of the seismic-trace data for the repetitions as disclosed in Weinstein et al. U.S. Pat. No. 3,946,357 by means of true amplitude summation would be dominated by the noise.
Another relatively simple approximation results when the assumption is made that s.sub.i is the same for all seismic-trace signals and that n.sub.i.sup.2 is approximately equal to the average absolute value of the i-th trace. Then, Equation (12) reduces to: ##EQU9## This process is comparable to applying "digital AVC" to each trace before vertical stacking and is related to mantissa-only, sign bit, and automatic gain control (AGC) recording implemented in various commercially available digital field recorders.
However, again consider several repetitions of a swept frequency "chirp" signal at each of a plurality of vibration points for providing a set of seismic-trace data from which the true earth response is to be statistically estimated. Assume one recording within the set to be in the presence of very high burst noise. Assume also that the very noisy recording can be multiplied by a weight between one and zero. If, on the one hand, the weight is one, the result is true amplitude summation and, in such a case, the one noisy recording will overwhelm the others, the resulting estimate being that of noise only as indicated above. In contrast, a factor less than one can be applied for weighting the very noisy recording so that the impact on the estimate is approximately the same as the less noisy recordings, which is the motivation for mantissa-only, sign bit, and AGC recording. Importantly, such features can be implemented on almost all field hardware inexpensively in a way which is time and trace variable (which is critical since noise is time and trace variable). However, mantissa-only and sign bit recording affect the seismic-trace data frequency content due to stepwise transitions in the digitized traces and thereby cause a loss of informational content. AGC recording does not generally result in a loss of informational content, but at best AGC recording reduces the impact of the very noisy recording on the estimate only slightly, since the noise within the recording will dominate that recording as there is no noise suppression within the recording, that is, the recording is virtually all noise. Hence, the noise is significant because the very noisy recording is accorded the same contribution to the estimate as the less noisy recordings.
Embree U.S. Pat. No. 3,398,369 discloses an approximation known as "diversity stacking" which assumes that the true earth response seismic signal amplitude can be estimated from the total powder in the "early" portion of the trace (that is, s.sub.i .alpha..sqroot..SIGMA.X.sub.i,j.sup.2) and the noise power from the total power in the "late" portion of the trace (that is, n.sub.i.sup.2 .alpha..SIGMA.X.sub.i,j.sup.2). Equation (12) can then be estimated from the ratio of these true earth response and noise power estimations under the assumption that s.sub.i is nearly the same for all seismic-trace signals and variations in the total power from one trace to the next for any given window are dominated by variations in noise. Accordingly, Equation (12) can be approximated by: ##EQU10## Implementation of diversity stacking consists of first partitioning each seismic-trace signal into a series of windows. Next, the total power in each window is computed, and the seismic-trace data is scaled by a windowwise linear function of the inverse of the power in that window and the power in the previous window. (It should be noted, in passing, that the calculation and application of weighting scale factors could also be accomplished over "moving windows;" such a scheme, however, would require more computational complexity.) Finally, the scaled seismic-trace data is algebraically summed and normalized prior to recording. The normalization scale factors are inversely proportional to the sum of the weighting scale factors on a per window basis.
Diversity stacking is time variable depending on the length of the portion of the recording on which the weighting scale factor is based. This process is implemented in various commercially available digital field recorders and has been used for reducing burst and spike noise in vibroseis seismic prospecting recordings.
Nevertheless, in some commercially available digital field recorders wherein diversity stacking has been implemented, the weighting scale factors are determined by: ##EQU11## where c.sub.1, c.sub.2, c.sub.3, .epsilon..sub.1, and .epsilon..sub.2 are constants and comparative limits. This reduces the chance that a near-zero seismic-trace signal will dominate the vertically stacked recording and at the same time allows a high-noise trace to be "muted" or "blanked."
However, diversity stacking disclosed in Embree U.S. Pat. No. 3,398,396 is based on the assumption that the magnitude of the noise on each recording is calculable so that the exact weighting scale factor can be calculated for each recording at a given vibration point that will maximize the signal-to-noise ratio. (The weighting scale factor is the inverse square of the noise amplitude.) Unfortunately, one does not know the amplitude of the noise on each recording. One only knows the amplitude of the noise plus the true earth response seismic signal.
In any event, diversity stacking has been found to have various limitations. One limitation is that the use of diversity stacking does not yield an optimum signal-to-noise ratio in circumstances where the level of the true earth response seismic signal is comparable in magnitude to the nonstationary noise. Since the true earth response seismic signal, which includes components such as ground roll and refractions, is often comparable in magnitude to burst noise which appears during vibroseis seismic prospecting, diversity stacking does not always provide adequate noise suppression in noisy vibroseis seismic prospecting. Another limitation of diversity stacking is the difficulty and complexity of implementation in field hardware.
Reject recording is implemented in some commercially available digital field recorders. Reject recording causes a weighting scale factor of zero to be applied if a predetermined threshold is reached. The effect of noisy recordings are eliminated, that is, rejected. However, reject recording affects the seismic-trace data frequency content and thereby causes a loss of informational content, and furthermore, the required adjustments of threshold for producing such rejection have proven difficult to carry out in the field. Reject recording is dependent upon predetermined selection of threshold, which, if improperly selected, can, on the one hand, completely eliminate all true earth response seismic signals or, on the other hand, fail to reject any noise.
Other noise suppression schemes have been proposed which are independent of weighting prior to vertical stacking. Examples of such schemes are disclosed in Schmitt U.S. Pat. No. 3,744,019 and Siems U.S. Pat. No. 3,894,222, for instance.